Modelo de Regresión Logística

Debido a que ya tenemos los datos por las entregas anteriores procedemos a elaborar un modelo de regresión logística para la variable “EsCara” utilizando validación cruzada.

# install.packages("caret")
# install.packages("e1071")  # requerido por caret para modelos SVM y otros

library(caret)
library(e1071)

## Warning: package 'e1071' was built under R version 4.4.3

set.seed(randomSeed)  # Para reproducibilidad

control <- trainControl(method = "cv",    # cross-validation
                        number = 10,      # número de folds
                        classProbs = TRUE,  # para clasificación
                        summaryFunction = twoClassSummary)  # para métricas como ROC

# Puedes cambiar a Accuracy, Sensitivity, etc.
#modeloBarata <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
#modeloMediana <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
modeloCara <- train(EsCara ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
                   data = train_data,
                   method = "glm",
                   family = "binomial",
                   trControl = control,
                   metric = "Accuracy")

## Warning in train.default(x, y, weights = w, ...): The metric "Accuracy" was not
## in the result set. ROC will be used instead.

#print("Modelo de EsBarata")
#print(modeloBarata)

#print("Modelo de EsMediana")
#print(modeloMediana)

print("Modelo de EsCara")

## [1] "Modelo de EsCara"

summary(modeloCara)

## 
## Call:
## NULL
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -8.4410     1.0141  -8.323  < 2e-16 ***
## LotArea       0.2312     0.2635   0.878 0.380169    
## OverallQual   2.0755     0.4132   5.023 5.09e-07 ***
## YearBuilt     0.6113     0.4224   1.447 0.147831    
## GarageCars    1.6028     0.4714   3.400 0.000673 ***
## GrLivArea     1.0173     0.3026   3.362 0.000775 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 369.23  on 1023  degrees of freedom
## Residual deviance: 129.42  on 1018  degrees of freedom
## AIC: 141.42
## 
## Number of Fisher Scoring iterations: 9

Analizando los coeficientes podemos ver que varias de las variables utilizadas tienen un p-value menor a 0.05, pero no son todas, LotArea y YearBuilt tienen un p value demasiado elevado, lo que me lleva a pensar que realmente no necesariamente se correlacionan con el precio de la vivienda.

Utilizando el modelo con el conjunto de verificación podemos ver que:

predicciones <- predict(modeloCara, newdata = test_data)

# Convertir a clases (0 o 1) usando un umbral (por ejemplo, 0.5)
# predicciones <- ifelse(probabilidades > 0.5, 1, 0)
confusionMatrix(as.factor(predicciones), as.factor(test_data$EsCara))

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  No  Sí
##         No 414   8
##         Sí   4  10
##                                           
##                Accuracy : 0.9725          
##                  95% CI : (0.9524, 0.9857)
##     No Information Rate : 0.9587          
##     P-Value [Acc > NIR] : 0.08708         
##                                           
##                   Kappa : 0.6109          
##                                           
##  Mcnemar's Test P-Value : 0.38648         
##                                           
##             Sensitivity : 0.9904          
##             Specificity : 0.5556          
##          Pos Pred Value : 0.9810          
##          Neg Pred Value : 0.7143          
##              Prevalence : 0.9587          
##          Detection Rate : 0.9495          
##    Detection Prevalence : 0.9679          
##       Balanced Accuracy : 0.7730          
##                                           
##        'Positive' Class : No              
##

Viendo los resultados del modelo, podemos ver que aunque el “Accuracy” es muy alto, realmente el modelo es medio malo, ya que nuestro Balanced Accuracy apenas llega a 77%, esto se debe a que realmente la cantidad de casas que cumplen nuestra definición de “cara” es extremadamente alta (2 desviaciónes estándar por encima de la media). Por lo tanto tenemos que balancear la muestra para que el modelo pueda aprender características sobre este conjunto de datos reducido.

# install.packages("ROSE")
library(ROSE)

## Warning: package 'ROSE' was built under R version 4.4.3

## Loaded ROSE 0.0-4

# Cargar datos
training_data <- read.csv("data/train.csv")

# Reemplazar valores NA con 0
training_data[is.na(training_data)] <- 0

mean_price <- mean(training_data$SalePrice, na.rm = TRUE)
sd_price <- sd(training_data$SalePrice, na.rm = TRUE)

# Definir los límites
lower_limit <- mean_price - (1 * sd_price)
upper_limit <- mean_price + (2 * sd_price)


# Crear la variable de clasificación
training_data$Category <- ifelse(training_data$SalePrice < lower_limit, "Baratas",
                           ifelse(training_data$SalePrice > upper_limit, "Caras", "Medianas"))
training_data <- training_data %>%
  mutate(
    EsBarata = ifelse(SalePrice < lower_limit, "Sí", "No"),
    EsCara = ifelse(SalePrice > upper_limit, "Sí", "No"),
    EsMediana = ifelse(SalePrice >= lower_limit & SalePrice <= upper_limit, "Sí", "No")
  )
training_data$Category <- as.factor(training_data$Category)
balanced_data <- ROSE(EsCara ~   OverallQual  + GarageCars + GrLivArea + Category, data = training_data)$data

predictors <- balanced_data %>% select(-Category)
response <- balanced_data$Category

# Normalización (Estandarización Z-score)
preProc <- preProcess(predictors, method = c("center", "scale"))
predictors_scaled <- predict(preProc, predictors)

# Semilla para reproducibilidad
set.seed(randomSeed) 

# Separar datos en entrenamiento (70%) y verificación (30%)
train_indices <- createDataPartition(response, p = 0.7, list = FALSE)

train_data <- predictors_scaled[train_indices, ]
test_data <- predictors_scaled[-train_indices, ]

library(caret)
library(e1071)

set.seed(randomSeed)  # Para reproducibilidad

control <- trainControl(method = "cv",    # cross-validation
                        number = 10,      # número de folds
                        classProbs = TRUE,  # para clasificación
                        summaryFunction = twoClassSummary)  # para métricas como ROC

# Puedes cambiar a Accuracy, Sensitivity, etc.
#modeloBarata <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
#modeloMediana <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
modeloMejoradoCara <- train(EsCara ~ OverallQual  + GarageCars + GrLivArea,
                   data = train_data,
                   method = "glm",
                   family = "binomial",
                   trControl = control,
                   metric = "Precision-Recall AUC")

## Warning in train.default(x, y, weights = w, ...): The metric "Precision-Recall
## AUC" was not in the result set. ROC will be used instead.

#print("Modelo de EsBarata")
#print(modeloBarata)

#print("Modelo de EsMediana")
#print(modeloMediana)

print("Modelo de EsCara")

## [1] "Modelo de EsCara"

predicciones <- predict(modeloMejoradoCara, newdata = test_data)

# Convertir a clases (0 o 1) usando un umbral (por ejemplo, 0.5)
# predicciones <- ifelse(probabilidades > 0.5, 1, 0)
confusionMatrix(as.factor(predicciones), as.factor(test_data$EsCara))

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  No  Sí
##         No 190  10
##         Sí  22 215
##                                           
##                Accuracy : 0.9268          
##                  95% CI : (0.8982, 0.9494)
##     No Information Rate : 0.5149          
##     P-Value [Acc > NIR] : < 2e-16         
##                                           
##                   Kappa : 0.8532          
##                                           
##  Mcnemar's Test P-Value : 0.05183         
##                                           
##             Sensitivity : 0.8962          
##             Specificity : 0.9556          
##          Pos Pred Value : 0.9500          
##          Neg Pred Value : 0.9072          
##              Prevalence : 0.4851          
##          Detection Rate : 0.4348          
##    Detection Prevalence : 0.4577          
##       Balanced Accuracy : 0.9259          
##                                           
##        'Positive' Class : No              
##

¡Podemos ver que este modelo se comporta de una manera mucho mejor al anterior! Aunque el accuracy normal disminuyó considerablemente el accuracy balanceado aumentó a 92%! Esto se debe principalmente a que nuestro modelo ya es capaz de identificar muchas más casas que sí son consideradas caras, pueso que las incluimos más seguido dentro del dataset.

Análisis de Overfitting/Underfitting

¡Para analizar el overfitting/underfitting del modelo necesitamos evaluarlo con los datos de entrenamiento y comparar sus resultados con respecto a los datos de verificación!

library(plotly)

## Warning: package 'plotly' was built under R version 4.4.3

## 
## Adjuntando el paquete: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

print("Usando data de entrenamiento")

## [1] "Usando data de entrenamiento"

predicciones <- predict(modeloMejoradoCara, newdata = train_data)

# Convertir a clases (0 o 1) usando un umbral (por ejemplo, 0.5)
# predicciones <- ifelse(probabilidades > 0.5, 1, 0)
# confusionMatrix(as.factor(predicciones), as.factor(train_data$EsCara))


curva <- learning_curve_dat(dat = train_data,
                            outcome = "EsCara",
                            proportion = seq(0.1, 1.0, by = 0.1),
                            test_prop = 0.3,
                            method = "glm",
                            metric = "Accuracy",
                            family = "binomial")

## Training for 10% (n = 71)

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Training for 20% (n = 143)
## Training for 30% (n = 215)
## Training for 40% (n = 286)
## Training for 50% (n = 358)
## Training for 60% (n = 430)
## Training for 70% (n = 501)
## Training for 80% (n = 573)
## Training for 90% (n = 645)
## Training for 100% (n = 717)

# Graficar la curva de aprendizaje
ggplotly(
  ggplot(curva, aes(x = Training_Size, y = Accuracy, color = Data)) +
  geom_smooth(se = FALSE) +
  labs(title = "Curva de Aprendizaje - Regresión Logística",
       y = "Accuracy", x = "Tamaño del conjunto de entrenamiento")
)

## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Como podemos ver, las curvas tanto de validación como de entrenamiento se siguen muy cercanamente en los rangos cercanos a 300 datos y de 500 en adelante. En el final, aunque no convergen se puede ver que sí se encuentran muy cercanos entre sí por lo tanto no hay Overfitting. Tampoco creemos que haya underfitting, puesto que aunque sí están juntas la mayoría del tiempo, el valor de accuracy es demasiado alto (mayor a 94%). Por lo que consideramos que el modelo realmente sí aprendió de forma correcta luego de aplicarle un resampling a los datos de entrada para que no tuviera despreciara minorías.

Eficiencia del modelo

Es importante determinar como podemos mejorar nuestro modelo y esto se puede hacer modificando los hiperparametros usados. En la implementación actual del modelo de regresión logística utilizando el método “glm” dentro de la función train de la librería caret, no se puede realizar ajustes automáticos de hiperparámetros. Esto se debe a que la función glm en R, implementa la regresión logísitca estándar, no posee hiperparámetros intrínsecos que puedan ser optimizados.

En este caso el único “tunning” que se podría hacer es volver a realizar la ingeniería de caractrísticas y cambiar la selección de las variables predictoras pero esto involucraría cambiar el resto de los modelos que se hicieron en entregas anteriores por lo que no se procedera con este método y se dejara el tunning estandar que ofrece el modelo.

# install.packages("ROSE")
library(ROSE)

# Cargar datos
training_data <- read.csv("data/train.csv")

# Reemplazar valores NA con 0
training_data[is.na(training_data)] <- 0

mean_price <- mean(training_data$SalePrice, na.rm = TRUE)
sd_price <- sd(training_data$SalePrice, na.rm = TRUE)

# Definir los límites
lower_limit <- mean_price - (1 * sd_price)
upper_limit <- mean_price + (2 * sd_price)


# Crear la variable de clasificación
training_data$Category <- ifelse(training_data$SalePrice < lower_limit, "Baratas",
                           ifelse(training_data$SalePrice > upper_limit, "Caras", "Medianas"))
training_data <- training_data %>%
  mutate(
    EsBarata = ifelse(SalePrice < lower_limit, "Sí", "No"),
    EsCara = ifelse(SalePrice > upper_limit, "Sí", "No"),
    EsMediana = ifelse(SalePrice >= lower_limit & SalePrice <= upper_limit, "Sí", "No")
  )
training_data$Category <- as.factor(training_data$Category)
balanced_data <- ROSE(EsCara ~   OverallQual  + GarageCars + GrLivArea + Category, data = training_data)$data

predictors <- balanced_data %>% select(-Category)
response <- balanced_data$Category

# Normalización (Estandarización Z-score)
preProc <- preProcess(predictors, method = c("center", "scale"))
predictors_scaled <- predict(preProc, predictors)

# Semilla para reproducibilidad
set.seed(randomSeed) 

# Separar datos en entrenamiento (70%) y verificación (30%)
train_indices <- createDataPartition(response, p = 0.7, list = FALSE)

train_data <- predictors_scaled[train_indices, ]
test_data <- predictors_scaled[-train_indices, ]

library(caret)
library(e1071)

set.seed(randomSeed)  # Para reproducibilidad

control <- trainControl(method = "cv",    # cross-validation
                        number = 10,      # número de folds
                        classProbs = TRUE,  # para clasificación
                        summaryFunction = twoClassSummary)  # para métricas como ROC

# Puedes cambiar a Accuracy, Sensitivity, etc.
#modeloBarata <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
#modeloMediana <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
#                   data = train_data,
#                   method = "glm",
#                   family = "binomial",
#                   trControl = control,
#                   metric = "ROC")
modeloMejoradoCara <- train(EsCara ~ OverallQual  + GarageCars + GrLivArea,
                   data = train_data,
                   method = "glm",
                   family = "binomial",
                   trControl = control,
                   metric = "Precision-Recall AUC")

## Warning in train.default(x, y, weights = w, ...): The metric "Precision-Recall
## AUC" was not in the result set. ROC will be used instead.

#print("Modelo de EsBarata")
#print(modeloBarata)

#print("Modelo de EsMediana")
#print(modeloMediana)

print("Modelo de EsCara")

## [1] "Modelo de EsCara"

predicciones <- predict(modeloMejoradoCara, newdata = test_data)

# Convertir a clases (0 o 1) usando un umbral (por ejemplo, 0.5)
# predicciones <- ifelse(probabilidades > 0.5, 1, 0)
confusionMatrix(as.factor(predicciones), as.factor(test_data$EsCara))

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  No  Sí
##         No 208  22
##         Sí  18 188
##                                           
##                Accuracy : 0.9083          
##                  95% CI : (0.8772, 0.9336)
##     No Information Rate : 0.5183          
##     P-Value [Acc > NIR] : <2e-16          
##                                           
##                   Kappa : 0.8161          
##                                           
##  Mcnemar's Test P-Value : 0.6353          
##                                           
##             Sensitivity : 0.9204          
##             Specificity : 0.8952          
##          Pos Pred Value : 0.9043          
##          Neg Pred Value : 0.9126          
##              Prevalence : 0.5183          
##          Detection Rate : 0.4771          
##    Detection Prevalence : 0.5275          
##       Balanced Accuracy : 0.9078          
##                                           
##        'Positive' Class : No              
##

Viendo la matriz de confusión obtenida vemos que estamos obteniendo más errores al predecir que una casa es cara, teniendo 22 casos (falsos negativos) que las casas si eran caras pero las clasifico como “No”. Luego vemos que que hubo 10 casos (falsos positivos) que fueron mal clasificados. Entendemos que estos errores implican en que se afecta en la eficiencia de las ventas demostrando posibles falsos datos a los compradores haciendo que exista posibles pérdidas en las oportunidades de venta.

profvis({
  # Cargar datos
  training_data <- read.csv("data/train.csv")

  # Reemplazar valores NA con 0
  training_data[is.na(training_data)] <- 0

  mean_price <- mean(training_data$SalePrice, na.rm = TRUE)
  sd_price <- sd(training_data$SalePrice, na.rm = TRUE)

  # Definir los límites
  lower_limit <- mean_price - (1 * sd_price)
  upper_limit <- mean_price + (2 * sd_price)


  # Crear la variable de clasificación
  training_data <- training_data %>%
    mutate(
      EsBarata = ifelse(SalePrice < lower_limit, "Sí", "No"),
      EsCara = ifelse(SalePrice > upper_limit, "Sí", "No"),
      EsMediana = ifelse(SalePrice >= lower_limit & SalePrice <= upper_limit, "Sí", "No")
    )
  training_data$Category <- as.factor(ifelse(training_data$EsCara == "Sí", "Cara",
                                             ifelse(training_data$EsBarata == "Sí", "Barata", "Mediana")))
  balanced_data <- ROSE(EsCara ~ OverallQual + GarageCars + GrLivArea + Category, data = training_data)$data

  predictors <- balanced_data %>% select(-Category, -EsBarata, -EsMediana, -EsCara)
  response <- balanced_data$EsCara

  # Normalización (Estandarización Z-score)
  preProc <- preProcess(predictors, method = c("center", "scale"))
  predictors_scaled <- predict(preProc, predictors)

  # Semilla para reproducibilidad
  set.seed(randomSeed)

  # Separar datos en entrenamiento (70%) y verificación (30%)
  train_indices <- createDataPartition(response, p = 0.7, list = FALSE)

  train_data <- predictors_scaled[train_indices, ]
  test_data <- predictors_scaled[-train_indices, ]
  train_response <- response[train_indices]
  test_response <- response[-train_indices]

  set.seed(randomSeed)  # Para reproducibilidad

  control <- trainControl(method = "cv",    # cross-validation
                           number = 10,      # número de folds
                           classProbs = TRUE,  # para clasificación
                           summaryFunction = twoClassSummary)  # para métricas como ROC

  # Puedes cambiar a Accuracy, Sensitivity, etc.
  # modeloBarata <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
  #                       data = train_data,
  #                       method = "glm",
  #                       family = "binomial",
  #                       trControl = control,
  #                       metric = "ROC")
  # modeloMediana <- train(EsBarata ~ LotArea + OverallQual + YearBuilt + GarageCars + GrLivArea,
  #                       data = train_data,
  #                       method = "glm",
  #                       family = "binomial",
  #                       trControl = control,
  #                       metric = "ROC")
  modeloMejoradoCara <- train(EsCara ~ OverallQual + GarageCars + GrLivArea,
                              data = train_data,
                              method = "glm",
                              family = "binomial",
                              trControl = control,
                              metric = "ROC") # Cambié la métrica a ROC para twoClassSummary

  # print("Modelo de EsBarata")
  # print(modeloBarata)

  # print("Modelo de EsMediana")
  # print(modeloMediana)

  print("Modelo de EsCara")
  print(modeloMejoradoCara)

  predicciones_prob <- predict(modeloMejoradoCara, newdata = test_data, type = "prob")
  predicciones <- ifelse(predicciones_prob$Sí > 0.5, "Sí", "No")

  confusionMatrix(as.factor(predicciones), as.factor(test_response))
})

## profvis: code exited with error:

En lo que respecta al tiempo y la memoria consumida se ve que para este modelo estos vaores son significativamente bajos. Se utilizo la libreria profviz para analizar estos datos y se obtuvo que

ReporteArboles

Flavio Galán, José Prince

2025-04-12

Modelo de Regresión Logística

Análisis de Overfitting/Underfitting

Eficiencia del modelo